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One of the popular modifications to the theory of general relativity is non-dynamical Chern-Simons 
(CS) gravity, in which the metric is coupled to an externally prescribed scalar field. Setting accurate 
constraints to the parameters of the theory is important owing to their implications for the scalar 
field and/or the underlying fundamental theory. The current best constraints rely on measurements 
of the periastron precession rate in the double-binary-pulsar system and place a very tight bound on 
the characteristic CS lengthscale fc" 1 < 3 x 10~ 9 km. This paper considers several effects that were 
not accounted for when deriving this bound and lead to a substantial suppression of the predicted 
rate of periastron precession. It is shown, in particular, that the point mass approximation for 
extended test bodies does not apply in this case. The constraint to the characteristic CS lengthscale 
is revised to k^ 1 < 0.4 km, eight orders of magnitude weaker than what was previously found. 
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I. INTRODUCTION 

Einstein's theory of general relativity (GR) has so far 
passed all observational tests with flying colors (for a re- 
view, see for example Ref. [1]). It is expected, however, 
that GR is the low-energy limit of a more fundamental 
theory unifying all forces of nature. In that case it is 
likely that higher-order curvature corrections exist in the 
theory, the effect of which may become apparent in the 
strong field regime. Chern-Simons (CS) gravity J5[ |3] is 
an example of such a higher-order modification to GR, 
in which the metric is coupled to a scalar field through 
the parity-violating Pontryagin density. In dynamical 
CS gravity, the scalar field is itself coupled to the met- 
ric through a wave equation sourced by the Pontryagin 
density. In non-dynamical CS gravity, the subject of the 
present work, the scalar field is externally prescribed, and 
is typically assumed to be a homogeneous field that only 
depends on cosmic time. The implications of CS gravity 
have been investigated in several astrophysical and cos- 
mological scenarios. Ref. [I] considered its effects on the 
polarization of the cosmic microwave background, iden- 
tifying the scalar field with the inflaton. Ref. 2j worked 
out the linearized theory and the propagation of gravi- 
tational waves. Ref. [5j suggested that CS gravity may 
provide a mechanism for creating the observed cosmic 
matter-antimatter asymmetry. Ref. [5] investigated the 
post-Newtonian expansion of CS gravity in the point- 
particle limit. An important feature of the theory is that 
it leads to a change of frame dragging effects around ro- 
tating objects, which can be used to constrain CS gravity 
[H[7]. Smith et al. [8], hereafter SE08, have calculated 
the CS modification to the gravitomagnetic field around 
a non-relativistic, constant-density spinning body. Using 
the measurement of the Lense-Thirring drag around the 
Earth by the LAGEOS satellites [S], they have set the 
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first constraint on the characteristic CS lengthscale (de- 
fined below), fc^ 1 < 1000 km. More recently, Yunes and 
Spergel [TP] , hereafter YS09, used measurements of the 
rate of periastron precession in the double pulsar system 
PSR J0737-3039 A/B [HI H2 to place a stringent con- 
straint on the CS lengthscale, fc" 1 < 3 fim, eleven orders 
of magnitude stronger than the previous bound. 

One of the outcomes of the calculation of SE08 is that 
the CS modification to the gravitomagnetic field is oscil- 
lating in space with a wavelength 27r/c cs 1 . This oscillatory 
character is due to the higher-order nature of the theory. 
In this paper, we show that the oscillations of the grav- 
itomagnetic field lead to a large suppression of the peri- 
astron precession rate due to three effects that were not 
considered in previous studies: (i) The constant-density 
approximation is relatively accurate for the Earth but not 
for a neutron star, for which the surface density is 7 or- 
ders of magnitude lower than the mean density. We show 
that this leads to a suppression of the gravitomagnetic 
field outside the source by a factor ~ 5/(k cs Rj{), where 
Ra is the radius of the rapidly rotating star A. The reason 
of this suppression is that less sharp boundaries result in 
a smaller excitation of high-frequency spatial oscillation 
modes, (ii) The test body used for the constraint (star 
B) is extended, rather than a point particle. If the ra- 
dius of the star Rb is larger than the CS wavelength, this 
leads to a suppression of the average force per unit mass, 
by a factor ~ 15/(k cs R B ) 3 . (Hi) Even though the eccen- 
tricity of the system is e ~ 0.09 <§; 1 [T3], the semi-major 
axis a is large enough that the orbital separation varies 
over the scale of many CS wavelengths. We show that 
this leads to a suppression of the secular rate of perias- 
tron precession by a factor ~ l/^/k cs ae. Mostly because 
of effect (ii), we expect that one cannot constrain the CS 
lengthscale to be much smaller than the size of a neutron 
star with the system considered. Properly accounting 
for these three effects, we indeed derive the revised con- 
straint k^ 1 < 0.4 km, eight orders of magnitude weaker 
than what was previously found, opening up the space of 
allowed values for the parameters of non-dynamical CS 
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gravity. Maybe more importantly, the invalidity of the 
point-particle approximation for the test body shows that 
the "effacing principle" which holds in standard GR (see 
for example Ref. [13] ) is not necessarily valid in alternate 
theories of gravity. 

This paper is organized as follows: In Sec. [TTJ we review 
the basic equations of modified CS gravity and define our 
notation. We then describe the resulting gravitomagnetic 
field around a spinning neutron star in Sec. Ill (the de- 
tailed calculations can be found in Appendix | A[) . Sec. IV 
outlines how the anomalous periastron precession rate is 
obtained and the effect of previous assumptions (we also 
provide a derivation of the rate of change of orbital el- 
ements in Appendix [fj]). Sees. |v| and VI constitute the 
core of this paper, where we compute the suppression of 
the anomalous precession rate when correctly accounting 
for the finite extent of the test body and the eccentric- 
ity of the orbit, respectively. We derive the new con- 
straint in Sec. |VII| and discuss how to accurately do so 
when the magnitude of the predicted effect is oscillatory 
in Appendix[Cj We conclude and mention future research 
directions in Sec. IVIIII 



z?, V and £ 2 S appropriately). These two cases constitute 
two qualitatively different theories. 

Requesting the action to be stationary under varia- 
tions of the metric results in the modified Einstein field 
equation: 



<V + 16^ S C M „ = 8^ [T- at + ] 



(3) 



where is the Einstein tensor, T?J at is the matter 
stress-energy tensor, 

C>» = d a d e° a ^V a R v l + \v T {d„§) e a ^R% a (4) 

is a four-dimensional generalization of the Cotton- York 
tensor and 
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g^ a $V a $ - g^Vi-d) (5) 



is the stress-energy tensor for the scalar field 

If the scalar field is considered as dynamical, then vary- 
ing the action with respect to $ results in the equation 
of motion for the scalar field: 



II. CHERN-SIMONS GRAVITY: 
EQUATIONS 



BASIC 



For a review on Chern-Simons modified gravity, we re- 
fer the reader to Ref. [3] . Here we simply recall the main 
equations and define our notation. We use geometric 
units throughout the paper. 

We consider the following action defining the modified 
theory: 



-g <Fx+ C 



-g d x 



tfRR^-g d 4 x 



-g d x. 



(1) 



In Eq. , the first term is the standard Einstein-Hilbcrt 
action and the second term corresponds to the matter 
contribution with lagrangian density £ m at- The third 
term is the CS modification [2 , which depends on the 
constant t 2 s (with dimensions of length squared). It cou- 
ples the dimensionless CS scalar field ■& to the metric 
through the Pontryagin density RR, which is the con- 
traction of the Riemann tensor R a p^ v and its dual: 



RR = \e^ aT R af) R 



(2) 



where t^ vcrr is the four-dimensional Levi-Civita tensor. 
Finally, the last term of Eq. (JlJ contains the canoni- 
cal scalar field lagrangian, with potential V(§). In the 
non-dynamical version of the theory initially suggested 
in Ref. [2 , /3 = 0. In dynamical CS theory, /3 = 1 (a non- 
zero ft can always be set to unity provided one rescales 



pn<d 



RR + pV(d), 



(6) 



where □ is the usual covariant d'Alembertian operator, 
□ ee j^V^Vv. Equation Q can also be obtained by 
taking the divergence Of Eq. (§, using the fact that 
V^C^y = —^dv'&RR, as well as the Bianchi identity 
V^G^j, = and the conservation of the matter stress- 
energy tensor V M T^ at 



0. 



We see that in the dynamical theory (f3 — 1), Eq. Q is 
an evolution equation for the scalar field, which is sourced 
by the Pontryagin density. In the non-dynamical theory 
(J3 = 0), Eq. M is a constraint equation, RR = (the 
Pontryagin constraint), and the scalar field i? is uncon- 
strained and needs to be externally prescribed. The Pon- 
tryagin constraint imposes a severe restriction on the set 
of allowed spacetimes in non-dynamical CS gravity. 

In this work, we follow the approach of Refs. [HI HO], in 
that we technically work in the dynamical theory, with 
simplifying assumptions that allow us to make a specific 
choice for the scalar field, casting this work rather in 
the frame of non-dynamical CS gravity. We work in the 
weak field and slow rotation approximation. Outside a 
spinning source of angular momentum J, the Pontryagin 
density is [HI [16] RR = 288A/J ■ r/r 8 . If we work to 
first order in M/r and J/r 2 , the Pontryagin density is 
therefore formally a second order term. For the sake of 
simplicity, and following SE08 and YS09, we assume a 
vanishing potential, V(i!}) — 0. The evolution equation 
for the scalar field is therefore 



™ « 0, 



(7) 



to lowest order in M/r and J/r 2 , i.e. i9 is approxi- 
mately a harmonic function of the coordinates. Follow- 
ing Refs. [1 [HI [TO], we make the choice ■& = = M, 
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the assumption being that ■d is some cosmological field 
that traces the evolution of cosmic time. We finally note 
that the energy density of the scalar field has to be much 
smaller than the cosmological mean density [5] , and must 
therefore be much smaller than that of the source. We 
therefore have T*„ < T^ at in Eq. Even though we 
are formally working in the dynamical theory, the ap- 
proximations made and our choice of scalar field, identi- 
cal to the "canonical choice" suggested in Ref . [5] , there- 
fore rather cast this work in the frame of non-dynamical 
CS theory, while allowing a small violation of the Pon- 
tryagin constraint (see also the discussion in Ref. [3]). 

With these approximations the observables only de- 
pend on the combination t^J), which has dimensions of 
length. We define the characteristic CS wavenumber, 
which has units of inverse length: 



k cs ee (8n£ 2 Jy 



(8) 



The correspondence with the notation of Refs. [51 [TU] is 
k cs = m cs = 2/t cs . 

We emphasize that we have made several simplify- 
ing and somewhat arbitrary choices (following previous 
works), so that one should consider the theory studied as 
a toy model, aimed at gaining some insight into the more 
complex underlying fundamental physics. 



III. GRAVITOMAGNETIC FIELD AROUND A 
SPINNING OBJECT 

We follow SE08 and work in the non-relativistic, slow 
rotation regime. We consider the stationary (time- 
independent) problem. We work with the usual grav- 
itomagnetic vector potential A, in the Coulomb gauge 
V • A — [gauge freedom on spatial coordinates allows 
this choice; in the full time- dependent problem gauge 
freedom is limited by the choice of $ which determines 
the time coordinate]. With this setting, SE08 have shown 
that the only difference between CS gravity and GR is 
Ampere's law, which takes the form 



1+ fcr^v x A 



-4ttj, 



(9) 



or, equivalently, in terms of the gravitomagnetic field B 
Vxl, 



VxB - k-^B = Anj, 



(10) 



where j is the mass current of the source. The standard 
GR equation is recovered for k cs — > oo. For non-zero 
fc~ 1 , Eq. ([9J is a third order partial differential equation, 
and its solutions are therefore qualitatively different than 
those of GR. 

SE08 have solved Eq. ^ and computed the resulting 
gravitomagnetic field inside and outside a constant den- 
sity object. This is appropriate for the Earth, in which 
the density varies by a factor of a few from the center 



to the edge. This is however not accurately represent- 
ing a neutron star, where the density varies smoothly 
from the core to the edge, where it is ~ 10~ 7 smaller 
than the mean density. This smooth edge decreases 
the CS-induced gravitomagnetic field outside the star for 
large k cs R, as the high-spatial-frequency modes pick up 
a smaller amplitude. In Appendix |A"} we solve the mod- 
ified Ampere's equation for a more realistic neutron star 
density profile: 



p(r) = p c [1 - {r/Rf 



(11) 



We compute the gravitomagnetic field and obtain that 
the CS modification Bcs, = B — Bgr, where Bgr is the 
gravitomagnetic field for standard GR, is of the form [8J 

B cs = 47rpi? 2 j£>i(r)fi + D 2 (r)r x ft 

+ %)rx(fx(i)} J (12) 

where p = 3M/(4irR 3 ) is the mean density of the star. 
Outside the star, the functions Di are rescaled from the 
result of SE08 by a factor of 5fa(k cs R)/[k cs Rj>i(k cs R)] 
and are given by 



Di{r > R) = Wj 3 (k cs R) 



yi(k cs r) 



kcsT 

D 2 (r > R) = 5j 3 (k CB R) yi (k cs r), 
D 3 {r >R) = 5j 3 (k cs R)y 2 {k cs r), 



(13) 

(14) 
(15) 



where j n and y n are the order-n spherical Bessel func- 
tions of the first and second kind, respectively. We find 
that for k cs R 3> 1, the envelope of the oscillating grav- 
itomagnetic field outside the star is reduced by a factor 



5/(/e cs i?) when using the density profile (111 rather than 
assuming a constant density. We will show in Sec. |VII| 
that our final constraint is weakly dependent on the ex- 
act density profile and the simple form ( 11 ) is sufficiently 



accurate for the purpose of this calculation. 



IV. ANOMALOUS PERI ASTRO N PRECESSION 
IN CS GRAVITY 

In the binary pulsar, the orbital angular momentum 
is nearly aligned with the spin vector J of the central 
rotating source [T7], as one may expect from formation 
mechanisms. With this geometry, the CS gravitomag- 
netic field leads to an anomalous secular precession of 
the periastron with a rate (see Appendix [B] for a deriva- 
tion) : 



3(! _ e 2)3/2 



<r 3 J B (r-)[2+(e- 1 + e )cos/]) (16) 



where B (r) = A-KpR 2 [Di(r) — D 3 (r)], f is the true 
anomaly of the Keplerian trajectory, r(f) is the distance 
from the central rotating object and has the standard 
form given by Eq. (B13) for an ellipse of semi-major axis 
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a and eccentricity e, and (X)f is the angle- weighted av- 
erage of X over one orbit [see Eq. (B17)]. Moreover, we 
have used B (r), the mass- weighted average of the grav- 
itomagnetic field inside the test body (we will get back 
to this point shortly). 

In SE08 and YS09, the two following assumptions were 
implicitly made: 

(?) The test body was assumed to be a p oint mass, so 
that one can use B (r) = B Q (r) in Eq. (16 1. 

(ii) The eccentricity of the orbit was assumed to 
be small enough that the radius could be expanded 
around the semi-major axis inside the angular averag- 
ing, r w a(l — ecos/) and in particular Bq(t) ~ Bo(a) — 
aeB' (a) cos /. 

With these assumptions, one would obtain, to lowest 
order in eccentricity, cL>cs = 2[£?o( a ) — a ^o( a )l- 

Assumption (i) is only valid if the size of the test body 
is much smaller than the wavelength of the oscillating 
CS gravitomagnetic field, 27rfc~ 1 . This is indeed the case 
for the analysis done by SE08, who used the LAGEOS 
satellites [9] as test bodies, and who derived a relatively 
weak constraint fc" 1 < 1000 km. In this work as in YS09, 
the test body is star B of the double-binary pulsar, with 
radius Rb ~ 10 km. YS09 derive a constraint /c" 1 < 
3 /zm, many orders of magnitude smaller than the size 
of the test body. In that case, the CS gravitomagnetic 
force /cs = — 4u x Bqs oscillates a large number of times 
within the test body resulting in a near cancellation of 
the average force per unit mass. We illustrate this effect 
schematically in Fig. [T] We will show in Sec. [V] that the 
mass-averaged gravitomagnetic field in star B is 



Bo{r) = 15 7^F Bo(r) ' 



(17) 



which is suppressed by a factor ~ 15/(fc cs i?s) 3 from the 
point-mass case when fc cs i?s 3> 1. 

Let us now examine assumption [ii). The orbital sep- 
aration of the binary varies within the range a — ae < 
r < a + ae, i.e. the difference between the radial sep- 
aration at apoastron and periastron is 2ae. This cor- 
responds to 2aefc cs /(27r) oscillations of the CS gravito- 
magnetic field. In the case studied by SE08, a ~ 12000 
km, e < 0.01 and therefore 2ae <C 27r/c~ 1 when the con- 
straint of SE08 is saturated, i.e. the test body remains 
on roughly the same wavefront of the gravitomagnetic 
field throughout its orbit. In the double-binary-pulsar 
system, however, a s» 4 x 10 5 km and e « 0.09 [T2] . 
and therefore 2ae » 7 x 10 4 km » fc^ 1 , since SE08 have 
shown that fc^ 1 < 1000 km. The test body therefore goes 
through many oscillations of the gravitomagnetic field 
during each orbit, and the resulting force nearly averages 
out. We illustrate this effect in Fig. [2j Mathematically 
speaking, since B (r) oc cos(fc cs r), one cannot approxi- 
mate Bo(r) s» Bo(a) — aeB' (a) cos / unless k cs ae <C 1. 
We will show in Sec. |VI| that properly averaging the grav- 
itomagnetic force on an orbit leads to an additional sup- 
pression of the secular rate of periastron precession by a 
factor ~ 1/ (eV 'irk cs ae) . 



, |.| , 
h 



Star A 




FIG. 1. Schematic representation (not to scale) of the first ef- 
fect discussed in this paper. The solid and dotted lines repre- 
sent the crests and troughs of the oscillating gravitomagnetic 
field, with wavelength 27rfc~ 1 . The test body, star B, extends 
over many wavelengths of the CS gravitomagnetic field. The 
force per unit mass is therefore largely suppressed compared 
to the case of a point particle. 




FIG. 2. Schematic representation (not to scale) of the sec- 
ond effect discussed in this paper. The solid and dotted lines 
represent the crests and troughs of the oscillating gravitomag- 
netic field, with wavelength 2-7rfc~ 1 . The figure represents the 
orbit as seen face-on. During its eccentric trajectory, star B 
goes through many wavelengths of the oscillating CS gravito- 
magnetic field generated by the rotating star A. The net force 
nearly averages out, except for a small section at the pericen- 
ter and apocenter. Note that the scale used here is different 
from that of Fig. [l] 



Before turning to the computation of the effect men- 
tioned above, we emphasize that these two effects are ad- 
ditive, i.e. the mass-averaged gravitomagnetic field B$(f) 
remains a rapidly oscillating function of the position of 
the test body's center of mass [it is just uniformly sup- 
pressed by a constant factor as can be seen in Eq. (17)]. 



V. AVERAGE FORCE PER UNIT MASS FOR 
AN EXTENDED BODY 

In standard GR, all metric components are slowly vary- 
ing functions of position in the far-field of a given source. 
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As a consequence, the motion of an extended test body of 
characteristic size R at a distance r> JJ from a central 
source can be approximated by the motion of its center of 
mass, regardless of its internal structure. The latter only 
appears in corrections of relative amplitude (R/r) 2 . This 
"effacement" of internal structure [14], however, does not 
necessarily apply in alternative theories of gravitation 
[18j . in particular in non-dynamical CS gravity, as we 
shall show in this section. 

Consider a spherically symmetric test body of radius 
R, mass M and mass density p. The spherical Bessel 
functions can be written in terms of simple trigonometric 
functions, and the local CS gravitomagnetic force fcs = 
—4v x Bqs is °f the form 



fcs(r) = f c {r) cos(fc cs r) + f s (r) sin(& cs r), 



(18) 



where f c (r) and / s (r) are slowly varying functions of r 
(they are or the form 1/r"). Denoting r*o the position of 
the test body's center, the average force per unit mass 
on the entire body is 



/csM 



p(r')f(\r + f'\)d 3 r' 



1 

M 

« jjfM J p(r') cos (fc C8 |f + r'\) dV 

+ jjj, /.(ro) J P(r') sin (k cs \r + f'|) dV(19) 

where V is the body's volume, and since the amplitudes 
f c and f s remain nearly constant over the body's extent 
(if R <C r ), we can take them out of the integrals. We 
define the projected linear density 



A(r') 



P (Vr /2 +y 2 ) 2-Kydy. (20) 



Neglecting the curvature of planes of constant phase 
(valid if fc cs a -C a 2 /R 2 , indeed satisfied when our con- 
straint is saturated), we rewrite the average force per 
unit mass as 



1 



R 



/csM = Jjfc( r o) J A(/)cos[fc cs (r + r')]dr' 

i r R 

-J s (r a ) J ^X(r') S m[k CB (r + r')]dr'(21) 



M" 



Expanding the trigonometric functions and using the fact 
that A is an even function, we arrive at 



/csM =Z(k cs R)fcs(ro), 
where we have defined 



X(k cs R) 



M 



A(r') cos ( fc cs r')dr' . 



(22) 



(23) 



If the body's extent is much less than a CS wavelength, 
k cs R <§C 1, and X(k cs R) — 1 so 



/cs( r o; k cs R < 1) = /cs(»"o)- 



(24) 



This assumption was implicitly made in Refs. [5J QJj] to 
compute the force acting on the test body. 

For a general value of k cs R, however, I(k cs R) ^ 1, 
and in particular, for k cs R 3> 1, the rapidly oscillating 
integrand in Eq. p3| leads to a large suppression of X. 
We start by rewriting I, after integrating by parts and 
noticing that A(i?) = and d\/dr' = -2-nr' p(r'): 



T{k cs R) 



1 AttR 



k cs R M J 



r' p(r') s\a(k cs r')dr' 



k cs R Jq 



xp(x) sm(k cs Rx)dx, (25) 



where we have defined the normalized density profile 
p(r') = p(r') /p and made the change of variable x = 
r' I R. For an approximate neutron star density profile as 
in Eq. (Ill, we obtain 



I(fc cs i?) = 15 



h{k ca R) 
(k cs R) 2 ■ 



(26) 



We see that for k cs R » 1, X ~ 15/(fc cs i?) 3 and the aver- 
age force per unit mass is largely suppressed. Note that 
even for a constant density profile (and therefore sharp 
edges), one would obtain X(k cs R) = 3ji(k cs R) / (k cs R) ~ 
3/(fc cs i?) 2 for k cs R ^> 1, i.e. the main cause of the sup- 
pression is the fact that the test body is extended, and 
the smooth boundary leads to an additional suppression 
by a factor ~ b/(k cs R), just as in the case of the rotating 
neutron star that sources the CS gravitomagnetic field. 

We therefore conclude that the average force per unit 
mass on neutron star B, with radius Rb and situated at 
a distance r from the rotating body (star A) is 



/csM = 15 



j2(k cs R B ) 
(k cs R B ) 2 



/csO)- 



(27) 



This result can equivalently be quoted in terms of 
the mass-averaged gravitomagnetic field, Eq. (17 1, since 
fcs — — 4t7 x Bcs and v is a slowly varying function of 
position. 



VI. EFFECT OF ECCENTRICITY ON THE 
SECULAR CHANGE OF ORBITAL ELEMENTS 

Since we already know from Solar System constraints 
[8J that k cs a > 400 ^> 1 for the system considered, we can 
first expand the spherical bessel functions in k cs r ^> 1. 



Using Eqs. (17) and ( 13 )-( 15 ), with 



B (r) = 3 ^n A [D 1 (r)-D 3 (r)}, (28) 



we arrive at 



B (r) w B 



cos(fc cs r) 



(29) 



(> 



where we have defined the amplitude 



B a = -225-—il A ] 3 (k cs R A )— — — — 

Ra {KcsRBr 



(30) 



We now have all elements at hand to compute the rate 
of periastron precession. Using Eq. (29) in Eq. (16), we 
obtain, to lowest order in eccentricity: 



^cs 



4ff 



(cos[fc cs r(/)] cos /)/. 



(31) 



If it were not for the rapidly oscillating function, this 
term would in fact be independent of eccentricity to low- 
est order in e; however, we will see that for k cs ae 3> 1, 
the 1/e dependence remains so this is the dominant term 
in Eq. ((16). 



In the small eccentricity regime e <C 1, we have r(f) 
o(l — ecos/) and therefore 



(cos[fc cs r(/)]cos/)/ 



sin(fc cs a) (sin(fc cs ae cos /) cos /) / 
sin(fc cs a) Ji(fc cs ae), (32) 



where J\ is the order- 1 Bessel function of the first kind. 
In Eq. ( (32( ) we have expanded cos[r(/)] and used the 
fact that the term proportional to cos(fc cs ae cos /) cos / 
averages to zero (it flips sign under the change / — > it — 



/). Strictly speaking, Eq. (32) is only valid if k cs ae 2 <C 



1, required to approximate cos[fc cs r(/)] rs cos[fc cs a(l — 
ecos/)]; the final result is however independent of the 
order with which we take the limits k cs a 3> 1, e <C 1. 
Now expanding the Bessel function for k cs ae ^> 1, we 
finally obtain 



4B . ,, V2cos(fc cs ae+ f) 

wes » sm(fc cs a) . 

k cs ae \fiik cs ae 



(33) 



We see that the rate of periastron precession is sup- 
pressed by a factor ~ 1/ 1 \f"nk cs ae for k cs ae 3> 1. 

We can understand the magnitude of the suppression 
with the following simple argument. When k cs ae 1, 
the test body goes through a large number of peaks and 
troughs of Bqs during an orbit, the effect of which av- 
erages out nearly exactly. The only locations where the 
test body remains on the same wavefront of the oscil- 
latory metric is at the pericenter and apocenter of the 
trajectory, where the orbit osculates a circle (see Fig. [5]). 
We expand r(f) « a(l—e cos /) near the pericenter (r_ 
(V. 1 — e), /_ = 0) or apocenter (r + = a(l + e), / + = tt) 
and obtain 



r± ~ ^fae 



(/-/±) 2 



(34) 



The test body remains on the same wavefront of Bcs as 



long as k c 



r ± \ ^ 7r /2, which translates to 



|/-/±l 



< 



(35) 



Therefore the CS-induced forces do not average out only 
on a fraction of the orbit 



461 
2tt 



■\/7rfc cs ae' 



(36) 



and the secular variation of orbital parameters will be 
reduced by this factor with respect to a perfectly circular 
orbit. 



VII. REVISED CONSTRAINTS 

We can now use the measured periastron precession 
rate in the binary pulsar to set a constraint on the CS 
length fc^ 1 . Anticipating that the constraint will still be 
such that fc" 1 <C 10 km, we can expand Bo in k cs RA 3> 
1, k cs R,B 3> 1. Our final result for the rate of periastron 
precession is therefore 



u cs rj -900 

cos(k cs RA) sm(k cs R B ) 
R A {k cs R B f 



2 M A^ . n ,cos(fc cs oe+ \ 

-—n ASm (k cs a) {kcsae)3/2 



(37) 



This is an oscillating function of fc cs , which adds a lit- 
tle subtlety to the meaning of constraints (we discuss 
this point in Appendix [C]) . As a simple estimate, we can 
approximate wes by its envelope (i.e. replace all the si- 
nusoidal functions by unity): 



900 



2 M A 

7T R A 



n A (k cs ae)-^ 2 (k cs R B y 



(38) 



If one can measure the rate of periastron precession to 
be within an error 5uj from the value predicted by GR, 
then we can constrain the CS length to be less than 



kcs - UooV 2 Q A M A B[a ) 



1/5.5 



(39) 



With R B - 10 km, M A /R A ~ 1/5, o ~ 4 x 10 5 km, 
e ~ 0.09 and 2ir/Q, A w 22 ms, [H] and assuming (as in 
YS09) that the measured precession rate agrees with the 
GR prediction within the measurement error 8ui w 0.05 
degrees per year [TS], we obtain 



< 0.2 km 



5u) 

0.05 deg/yr 



1/5.5 



(40) 



This constraint is eight orders of magnitude weaker than 
what was found by YS09, but remains more than three 
orders of magnitude stronger than what can be inferred 
from Solar System tests. Note that because of the small 
exponent in Eq. (39), the constraint depends only weakly 



of the exact value chosen for Su> as well as on the exact 
value of the semi-major axis and eccentricity. Accounting 
more precisely for the oscillatory amplitude in Eq. (37) 
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leads to the formally more accurate constraint (see Ap- 
pendix [C] for details): 

fc^ 1 < 0.4 km [85% confidence level], (41) 



in the sense that there is still a 15% chance that fc^ 1 > 
0.4 km but the predicted u>cs remains consistent with 
observational constraints due to a low amplitude of the 
oscillatory functions in Eq. ( 37 1 . 



Finally, one might worry that the simple density profile 
that we have used for both neutron stars may not capture 
surface effects; however, for A: cs ~ 0.4 km, the size of the 
star is only a few times the CS wavelength and a coarse- 
grained density profile is largely sufficient. For example, 
using an even simpler constant density profile for both 
stars would only change the constraint by a factor of 
~ 2.5, making it tighter. We also explicitly evaluated the 
integral ( 25 ) for realistic density profiles computed with 
physically motivated equations of state (EOS) [19], and 
obtain 2" ~ 3/(fc cs i?) 2 5 for a wide range of neutron star 
masses and independently of the EOS chosen (the EOSs 
used where those of Refs.|20i and [2TI 22 ). We therefore 



conclude that using the density profile (11) provides a 
conservative, yet relatively accurate constraint, and we 



adopt Eq. (41) as our final result. 



VIII. CONCLUSIONS 

In this work, we have re-evaluated the correction to the 
rate of periastron precession in the double-binary-pulsar 
induced by the Chern-Simons modification to general rel- 
ativity. We found it to be largely suppressed compared to 
previous results, where several important effects were not 
accounted for. This large suppression is due to the os- 
cillatory character of the CS gravitomagnetic field which 
leads to near cancellations of the mean force inside an ex- 
tended test body and during an orbit, as well as a reduc- 
tion of the modified gravitomagnetic field generated by 
the source if its boundary is not sharp. We have revised 
the constraints on the CS characteristic length accord- 
ingly: fc" 1 = (8tt£^ s '&)~ 1 < 0.4 km. This new constraint 
opens up the space of allowed values for the coupling 
strength of the CS theory £l s or the scalar field deriva- 
tive One should keep in mind, however, that we have 
made several simplifying choices for the parameters of 
the theory (see Sec. [TT[ >, and this result only applies to 
the very specific simplified theory that was studied in this 
work. 

Finally, we would like to point out that a unique fea- 
ture of the CS gravitomagnetic field is the presence of a 
poloidal component, the effects of which have not been 
fully investigated yet. We show in Appendix |B 4| that 
it results in very interesting orbital dynamics. In par- 
ticular, it leads to a change of the angle between the 
orbital angular momentum and the spin of the source, as 
well as secular changes in the eccentricity and the magni- 
tude of the orbital angular momentum. Such effects are 
qualitatively different than standard GR predictions and 



could provide a powerful handle to constrain CS gravity 
further. A full study of the resulting dynamics would re- 
quire including simultaneously other spin-orbit coupling 
effects in GR and CS gravity [B], and will be the subject 
of future work. 
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Appendix A: Solution to the modified Ampere's 
equation in CS gravity 

In this Appendix we derive a solution to Ampere's 
equation in CS gravity, Eq. Such a solution was al- 
ready provided in SE08 for the case of a constant density 
body in solid rotation. We consider more general den- 
sity profiles, still assuming solid rotation with constant 
angular velocity fi. 

We start by writing the vector potential A = Aq R + 
Acs, where A GR is the solution of Ampere's law in stan- 
dard GR, 

V 2 I GR = -Anj. (Al) 
Equation ^ becomes, after multiplying by fc cs : 

V 2 [V x Acs + k cs A cs + V x A GR = 0, (A2) 

Requiring the term in brackets to vanish at infinity, the 
only solution for this Laplace's equation is 

V x Acs + k cs A cs = -V x 1 GK , (A3) 
or equivalently, 



B C s + k cs A C s = ~B GR , 



(A4) 



where Bqr — V x A GR an d B G s = V x Acs- We now 
work in spherical polar coordinates r, 9, <f>, choosing f2 
as the z-axis. The associated unit vectors are e r = f, 
eg = (sin0) _1 f x (f x fi) and = (sinfl) -1 ^ x f. We 
decompose a given vector V on this basis as V = V r e r + 
V e eg + V^e^. Since the GR gravitomagnetic field does 



not have any poloidal component (-Bqr = 0), projecting 
Eq. ( A4 ) on e$ gives us a first relation 



ut> — _u a<p 
- D cs ~~ ^cs^cS- 



(A5) 



Writing explicitly the curl in spherical polar coordinates, 
we also obtain (using azimuthal symmetry) 





1 


d 


B cs = 


r smO 






1 d 




B cs 


r dr 





(A6) 
(A7) 



We therefore see that the CS correction to the gravito- 
magnetic field is entirely determined by the poloidal com- 
ponent of the CS vector potential, A cs . We now multiply 
Eq. ( A3 ) by the operator (fc cs — V x ) . The Co ulom b gauge 
implies VxVxi = —V 2 A, and using Eq. ( Al ), we see 



that Acs is the solution of the inhomogeneous Helmholtz 
equation (similar to that obtained in SE08): 

V 2 l C s + k 2 cs A cs = 4ttj - fe cs B GR . (A8) 



We now project Eq. ( A8 1 on 
equation for A cs : 



and obtain the following 



+ 



d 2 


2 


d 


dr 2 


r 


dr 


1 




sin 


r 2 sin 


il 



AiTrp{r)VL sin 




(A9) 



where we used j = p(r)Q x f. This equation is separable, 



and the solution is of the form 



A cs (r,6) = iTTpR 3 ^ sin 8 A (^5«) 



(A10) 



where p = 3M/(4?ri? 3 ) is the mean density of the star. 
The function A(x) [with x = r/R] in Eq. (AlOl is the 



solution of the inhomogeneous spherical Bessel equation 



A"{x) + -A'{x) + 
x 



where we have defined 



A{x) = xp(x), (All) 



k cs R (A12) 



and p(x) = p(Rx)/p is the dimensionless density pro- 
file. Outside the star, the right-hand-side of Eq. (All) 
vanishes and A is a homogeneous solution: 



A(x > 1) = a ji(nx) + b yi(nx), 



(A13) 



where a and b are integration constants, and j n and y n 
are the order n spherical Bessel functions of the first and 
second kind, respectively. Expanding Eq. (All I near the 
origin, one finds that .4(0) = and A must be at least 



linear in x near x — 0. Denoting Ap a particular solution 
of Eq. (All I, we therefore have 

A(x < 1) =A P (x) +cji(kx), (A14) 

where we have eliminated the homogeneous solution pro- 
portional to yi as it diverges at x = 0. 

We therefore have three integration constants a, b, c to 
determine. Requiring A and its first derivative to be con- 
tinuous at the boundary of the star [as is required by the 
finiteness of all the coefficients in Eq. ( All |] provides only 
two constraints. The homogeneous solution (A13) is well 



behaved at infinity and one cannot, a priori, eliminate 
any of the constants a, b or combination thereof. 

This issue hints at potentially deeper problems with 
the theory, which has not yet been shown to be a well 
posed initial value problem (see however Ref. |23j which 
derived the boundary terms necessary to render the 
Dirichlct problem well defined). It is also possible that 
in a fully time-dependent study, the evolution selects a 
unique homogeneous solution outside the source. Since 
such a study is well beyond the scope of this work, and 
any choice of homogeneous solution outside the source 
does not affect any of our constraints quantitatively, we 
follow SE08 and choose the ansatz a = in Eq. (A13|. 



SE08 justify this choice by pointing out that in the limit 
k cs R <^ 1 and r > R, the function ji(k cs r) ~ ^k cs r in- 
creases with distance outside of the source, which they 
interpret as unphysical. 

With this ansatz, we can solve explicitly for the con- 
stants b and c given a particular solution Ap by requiring 
the continuity and smoothness of A at the boundary of 
the star. Solving the resulting linear second order al- 
gebraic equation, and using properties of the spherical 
Bessel functions, we arrive at 

b = K 2 j 2 ( K )A P {l) + Kj t (K) [A' P (l) - Ap(1)} , (A15) 

c - K 2 y 2 {K)Ap{\) + n yi {K) [A'p(l) - ^p(l)] .(A16) 

Before proceeding with specific examples, we write down 
the CS correction to the gravitomagnetic field in the same 
form as in SE08: 



4vrpi? 2 j£>i(r)0 + D 2 (r)f x fi 
+ D 3 (r)r x (f x 



are given by (with x 

D 1 (r)=2- 



where, using Eqs. (|A5j)-(A7) and (A10) 
rjR): 

Mx) 



D 2 (r) = kA{x), 
D 3 (r) = ^~ 



A'{x). 



In particular, outside the star, we obtain 
ing Eq. ( |A20| : 

2R 



(A17) 
the functions Di 

(A18) 
(A19) 
(A20) 
after simplify- 



Dx(r > R) 

D 2 {r>R) 
D 3 (r > R) 



b —yi(k cs r), 
r 

R y2{k cs r). 



b k, 

b k 



(A21) 

(A22) 
(A23) 
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For a constant density object, we have p(x) — 1 and 
the function Ap{x) — x/k 2 is a particular solution of 
Eq. (All ). We therefore obtain b — j2(«), c = f/2( K ) form 
Eqs. flA15p-dAl6l). One can then easily check that our 



solution for the gravitomagnetic field is then in agreement 
with that found by SE08. 

Assuming a constant density is a good approximation 
for the Earth, in which the density varies by only a fac- 
tor of a few from the center to the edge. In a neutron 
star, however, the surface density is in general of the or- 
der of ~ 10 times the mean density, and the density 
profile varies smoothly from the core to the edge. This is 
important as a smoother edge will result in a smaller am- 
plitude for the high spatial frequency modes k cs R S> 1, 
as we shall see below. 

As a simple yet more realistic model for a neutron star, 
we assume a density profile p(r) = p c [l — (r/R) 2 ] . The 
dimensionless density profile is then 



p(x) = -(l-x 2 ), 



(A24) 



where the normalization ensures that the average density 
is p. A particular solution for Eq. (All I with this density 
profile is 



A P (x) = 25- 



5i(l - a; 2 ) 
2 is 5 ' 



(A25) 



Using this solution into Eqs. (A15)-(A16) and after sim- 
plifications, we obtain 



6 = 5 



.73 («0 



.2/3(k) 



(A26) 



Using Eq. (A26) into Eqs. (A21>-(A23>, we obtain 



Eqs. ( 13 )-( 15 1 for the CS gravitomagnetic field outside 



a rotating neutron star. 



Appendix B: Derivation of the rate of change of 
orbital elements 

This appendix is dedicated to explicitly derive the rate 
of periastron precession given in Eq. ( 16 1 in the case 



where the orbital angular momentum is aligned with the 
spin vector of the central rotating object. We also show 
that the poloidal gravitomagnetic field leads to very in- 
teresting dynamics. We start with the general Gaussian 
perturbation equations. 



1. Gaussian perturbation equations 

The gaussian perturbation equations (see for example 
Ref. [24]) give the rate of change of osculating constants 
of motion q for a test mass evolving in a newtonian po- 
tential and subject to a perturbing force per unit mass 



dcj 
~dt 



<9q 
dv 



(Bl) 



Denoting p the total mass of the system, a trajectory in a 
newtonian potential is defined by 6 constants of motion: 
the energy per unit mass E — ^v 2 — p/r, the angular 

momentum per unit mass £ = r x v and the eccentricity 
vector e = -vx£ — r. Note that £ and e together contain 

only 5 constants of motion as £ ■ e = 0. The eccentricity 
vector has magnitude e, the eccentricity of the orbit, and 
points along the symmetry axis of the trajectory, towards 
the pericenter. Using Eq. (Bl ), we obtain that the rates 



of change of the constants of motion are 



dE 

~dt 

d£ 

dt 
de 

dt 



= v-6f, 
= fx 5f, 

= - \2(v-Sfy- (v-r)Sf- (r ■ 5f)v 
p L 



(B2) 

(B3) 
(B4) 



2. Application to a general gravitomagnetic force 

The gravitomagnetic force is 8f = —4v x B. This force 
does not do any work and therefore the energy (and as 
a consequence the semi-major axis a) is conserved, as 
can be seen from Eq. ( B2 ) . The gaussian perturbation 



equations for the angular momentum and eccentricity be- 
come: 



dt 
dt 

de _ 4 
dt p 



4(f • v)B - 4(f • B)v, 

(v-r)(vx B) + (B-£)v 



(B5) 
(B6) 



We now consider a general gravitomagnetic field of the 
form 

B = Bi(r)Q + B 2 (r)f xQ + B 3 (r)f x (f x Q), (B7) 



where Q is a fixed unit vector. Eq. (B7) describes 



both the GR gravitomagnetic field (with B±(r) = 
J/r 3 , £?2(r) = 0, By,{r) = |j/r 3 , where J is the spin 
of the rotating central object) and the CS gravitomag- 
netic field, Eq. (fl2|). 



We start by considering the simple case where f2 || £, 
which is the relevant case for the binary pulsar. We do so 
for a general field with no poloidal component (B2 = 0). 
We then separately consider the effect of the poloidal 
component for a general relative orientation of f2 and £. 



3. Case of Cl £, no poloidal component. 

In that case the gravitomagnetic fiel d sim plifie s to B = 
[-Bi(r) - B 3 (r)]Cl = B (r)l Eqs. |b5} and |B6} become 
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(using f ■ v = rr) 
dt 



ArrB Q {r% 
at 

de 4 r - 

= — Bo(r) rr(v x I) +£v 

fl L 



dt 



(B8) 
(B9) 



The right-hand-side of Eq. (B8) is a total time derivative 



and there are therefore no secular changes of the angular 
momentum vector: 



dt 

di/r 



= 0, 



(BIO) 



where 



(X) T = ^j T X(t)dt (Bll) 

is the time-average of X over one orbital period T . Since 
£ 2 = fj,a(l — e 2 ), no secular changes of I (and there- 
fore of its magnitude £) imply that there are no secular 
changes in the eccentricity e = ||e||. As a consequence, 
(de/dt) ■ e = [this can be also checked explicitly from 



Eq. (B9)]. Since the right-hand side of Eq. (B9) lies in 



the orbital plane, we conclude that the secular precession 
of the eccentricity vector must be of the form 



de\ 

di) T =UJU > 



(B12) 



where in the second line we have converted the time av- 
erage to an angular average 



1 



2tt 



(X) f = — / X(f)df. 
27T Jo 



(B17) 



We can apply this result to standard GR, with Bq(t) — 
— T^J/r 3 , and obtain the well known expression for the 
Lense-Thirring drag [note that this is the total rate of 
change of the longitude of the pericenter which is the 
sum of the rates of precession of the longitude of the 
ascending node and the argument of pericenter] : 



4J 



a 3(l_ e 2)3/2- 



(B18) 



For a slowly varying gravitomagnetic field, we can find 
the leading order of the precession rate in eccentricity, 
for e <C 1, by expanding r w a — ae cos /. We find 



u = 2[B (a)-aB' (a)] + O(e). (B19) 



4. Case of a poloidal gravitomagnetic field 

We now consider the case of a purely poloidal gravito- 
magnetic field, B = B{r)r x f2, for arbitrary orientation 
of the vector f2. This case is very interesting as it is com- 
pletely absent in standard G R an d, as we shall see, leads 
to unique effects. Using Eq. (B5), we find 



where u = £ x e. To get the rate of periastron precession 
uj, we need to project Eq. (B9) on the unit vector ii. 



Before doing so, it is useful to recall some properties of 
an elliptic orbit. The radial coordinate r is given by 



r(f) = Y 



P 



(B13) 



e cos / ' 

^ is the semi-latus rectum and / is the 



where p = a(l — e 
true anomaly, the angle between e and r 



Using conser- 



vation of angular momentum, I = r f, the derivative of 
r with respect to time can be written as 



r — —e sin /. 
P 



(B14) 



Finally, the projections of the velocity on the orthonor- 
mal vectors e, it are given by 



v ■ e 



1 ■ f 
--sin /, 

P 



-(e + cos/). 
P 



(B15) 



With these relations at hand and after simplifying, we 
finally obtain the rate of periastron precession: 



4(-B (r)(2 + e- 1 (l + e 2 )cos/)\ 



(1 - e 2 ) 3 / 2 



Mr) 2 



cos/ 



,(B16) 



dt 



= 4rrB(r)f x £1 



= 4r-esin/B(r) cos fie x ft) + sin/(w x .(B20) 
p L J 

The first term inside the brackets averages to zero (it is 
odd under the change / — > — /), and we therefore obtain 



1 / de \ ~ o 



(B21) 



7l = -(rB(r)sin 2 f) T . 



(B22) 



l\dt/T 

where we have defined the rate 

4 

-( 
P 

Note that 71 = 2B(a) to lowest order in eccentricity. 
Projecting this equation on £, we obtain that the mag- 
nitude of the orbital angular momentum changes with a 
secular rate 



(B23) 



This result is in stark contrast with GR, where (at least 
to lowest order) spin-orbit coupling does not change the 
magnitude of the orbital angular momentum or the spin, 
but only their relative orientation. We moreover obtain 
that the orientation of £ changes with the rate 



d£\ - - 

-} r = 7l e(^)e. 



(B24) 
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We can now readily obtain the rate of change of the ec- 
centricity, using £ 2 = (j,a(l — e 2 ): 



de\ 
~dt/T 



7l (l-e 2 )(fi-e). 



(B25) 



We see that circular orbits are unstable: a tiny eccentric- 
ity grows with rate 2B(a)(Cl ■ e). To obtain the rate of 
change of the orientation of the eccentricity vector, we 
first notice that e • I = at all times, which gives us 



(B26) 



Our last and most tedious task is to project Eq. (B6) on 
U = I x e. After some manipulations, we obtain 



de\ 

— ) ■ u 
dt/T 



72 



(ft.fi), 



(B27) 



where we have defined the rate 



72 = -(r J B(r)(2ecos/ + (l + e 2 )cos 2 /)) T , (B28) 
P 

which, to lowest order in eccentricity, has value 72 = 
25(a). Again, this indicates an instability: the eccen- 
tricity vector e — ee rotates with a divergent angular 
rate 2B(a)e~ 1 (Cl ■ u) for small eccentricities. 

We have therefore shown that the poloidal component 
of the gravitomagnetic field leads to unique dynamics: 

(i) A rotation of e and I around £ x e, with angular 
rate 7ie(f2 • £), vanishing for circular orbits. 

(ii) A precession of the eccentricity vector around £ 
with angular rate 72e _1 (f2 • u), divergent for arbitrarily 
small eccentricities. 

(Hi) A change of the magnitude of the orbital angu- 
lar momentum with rate £/£ = — 7ie(fi • e), accompa- 
nied with a change of the eccentricity with rate e = 
7i(l — e 2 )(f2 • e), non- vanishing even for arbitrarily small 
eccentricities. 

A full analysis of the dynamics of the system would re- 
quire accounting for the standard GR precession effects, 
spin-orbit coupling, etc., as well as a detailed analysis of 
the gravitomagnetic field induced by the orbital motion 
itself. This will be the subject of future work. 



Appendix C: Confidence level on the constraint 



Equation ( 37 1 is of the form 



lu cs = ak cs ■ sin(fc cs a) cos(fc cs ae + -) 

x cos(fc cs i?A) sin(fc cs i? s ), (CI) 

where a is a numerical constant depending on the system 
considered. Assuming Ra — Rb = R, we simplify this 
expression to 



— ' 
2' 



^CS = - fc cs 5 - 5 <S3, 



(C2) 



where S3 is the product of the 3 sinusoidal functions 

7T 

S3 = sin(fc cs a) cos(fc cs ae + — ) sin(2fc cs i?). (C3) 

From now on, to simplify the notation, we will use the 
CS characteristic length L = k^ 1 . L = corresponds to 
standard GR, and the CS-induced periastron precession 
rate is wes = f L 5 5 Ss. It will also be convenient to 
define the lengthscale Lq = (2Sdi/a) 1 ^ 5 ' 5 , so that one 
can rewrite 



L 



0.5 



8Cj. 



(C4) 



In Sec. VII we have, as a fir st approximation, simply 
used the envelope of Eq. (CI) |wcs| ~ aL 5 5 , to derive 
constraints on L from observational constraints |^csl ^5 
6lu. Clearly, however, there are always arbitrarily large 
values of L that will satisfy any observational constraint, 
for example L = Nna' 1 for any integer N, or values 
close enough from this. One should precisely quantify 
this caveat and rather quote constraints in the form U L < 
L max , with X% confidence", meaning that there is still 
a probability 1 — X% that |c!>cs| < Sui for L w L max . 

In order to do so, we need to evaluate the probabil- 
ity distribution for the magnitude of the product of the 
three sinusoidal functions S3. We first notice that the 
three functions can be assumed to be uncorrelated, be- 
cause of the largely different scales of their arguments 
(R -C ae <C a). Our first task is therefore to find the 
probability distribution for the amplitude of a single si- 
nusoidal function. Since we are mainly concerned about 
intervals where this amplitude might be small, we ap- 
proximate each sinusoidal function with a triangular pe- 
riodic function with the same tangents near zeros (see 
Fig. [3]). With this approximation, the probability that 
any single sinusoidal function has an absolute value in the 
range [e, e + de] is pi(e)de — ^de (this is correct to order 
0(e 2 ) for e <C 1 but wc formally use it for e € [0,n/2]). 
We now would like to compute the probability distribu- 
tion p 3 (e = I1S3I) for the absolute value of the product 
of the 3 uncorrelated sinusoidal functions. It is easier to 
work with itj = ln(2ej/7r), where q = | sin(...)| is the ab- 
solute value of one of the three sinusoidal functions. We 
have u — ln(2e/7r) = u\ + u-i + 113, where a change of 
variables gives us Pi(uj) = 2ei/7r = exp(itj). Therefore, 
we get 

Pa(u)= / Pi(ui)du! / p 1 (u 2 )du 2 pi{u - u x - u 2 ) 

= i-u 2 exp(u). (C5) 

We can now obtain the probability that the product of 
three sinusoidal functions is less than a certain value e: 



P(|5 3 | <e)= r P3 (e')de' 
Jo 



|("l + |ln(2 e /7r)| + i[ln(2 e /7r)] 2X 



(C6) 
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where we have used Eq. 




FIG. 3. To simplify the calculation, we approximate each 
sinusoidal function with a triangular periodic function with 
the same tangents near zeros. 
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(C5| for pz(u') and replaced 
Because of the (hie) 2 



u = ln(2e/-7r) in the final result 
term, we see that there is actually quite a significant 
probability that the product has a small amplitude. For 
example, there is a 48% chance that the product of sinu- 
soidal functions is less than 0.1 in magnitude and still a 
12% chance that it is less than 0.01. 

We can now precisely quantify how reliable a constraint 
is. For a given L, |wcs| < 5uj is equivalent to \S ?,\ < 
(Lq/L) 5 5 , the probability of which is given by Eq. (C6) 
with e = {L /L) 5 - 5 . We show the probability l-Pd^^ 
(Lq/L) 5 - 5 ) as a function of L in Fig. [ZJ We interpret 
this probability as a level of confidence on the constraint 



< L. 



For example, we obtain 

A;" 1 < 0.32 km 
fc" 1 < 0.54 km 



o confidence] 
o confidence! 



Even though here we have chosen the often used 68% and 
95% confidence intervals, we emphasize that the proba- 
bility distribution is not a gaussian. Moreover, we insist 
that the strict meaning of these confidence intervals is 
that one may still have fc" 1 equal to the maximum value 
quoted with the probability complementary to the confi- 
dence level. 



FIG. 4. Probability that the observational constraint |cjcs| < 
5u) is violated as a function of the CS lengthscale k^s ■ 
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